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Abstract — Wave propagation in ray-chaotic scenarios, charac- 
terized by exponential sensitivity to ray-launcliing conditions, is 
a topic of significant interest, with deep phenomenological impli- 
cations and important applications, ranging from optical com- 
ponents and devices to time-reversal focusing/sensing schemes. 
Against a background of available results that are largely focused 
on the time-harmonic regime, we deal here with short-pulsed 
wavepacket propagation in a ray-chaotic enclosure. For this 
regime, we propose a rigorous analytical framework based on 
a short-pulsed random-plane-wave statistical representation, and 
check its predictions against the results from finite-difference- 
time-domain numerical simulations. 

Index Terms — Ray chaos, short pulses, plane waves, random 
fields. 

I. Introduction and Background 

DURING the last decades, there has been a growing 
interest in the study of electromagnetic (EM) wave prop- 
agation in environments featuring ray-chaotic behavior, i.e., 
exponential separation of nearby-originating ray trajectories 
im, ||2|- Remarkably, the onset of such odd behavior is not 
necessarily a mere consequence of geometrical/structural com- 
plexity, but may occur even in deceptively simple coordinate- 
non-separable geometries with linear constitutive relationships 
as an effect of the inherent nonlinearity of the ray equation 
||3|. By exploiting the formal analogy with "billiards" and 
"pinballs" from classical mechanics (see, e.g., l^l-fTj for a 
sparse sampling, and [81 for a comprehensive review of the 
key concepts and ideas), simple examples of ray-chaotic closed 
(cavity-type) or open (scattering) propagation scenarios may 
be readily conceived (see, e.g., [9|-|17|). 

Besides the critical impUcations in the actual applicability 
of ray-based theory and related computational schemes (see, 
e.g., [10]), ray-chaotic manifestations have been shown to 
play a key role in several fields of application, including 
optical microcavities and lasers |[T8^-f2T1, characterization of 
scattering signatures [22J, ||23l . time-reversal-based focusing 
ll24l . Il25l and sensing ll26l . and reverberation chambers lIZTl - 
II291. 

From the phenomenological viewpoint, a large body of theo- 
retical, numerical and experimental studies have evidenced the 
presence of distinctive features (often with universal character) 
in the time-harmonic high-frequency wave dynamics observed 
in ray-chaotic scenarios, in terms, e.g., of the ensemble proper- 
ties of eigenvalues, eigenfunctions, and impedance/scattering 
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matrices, in close analogy with quantum chaos (see, e.g., ['301- 
[33] and references therein). These features are not observable 
in ray-regular (e.g., coordinate-separable) configurations, and 
are connected to the spectral properties exhibited by certain 
classes of random matrices |34|. Hence, although the (linear) 
wave dynamics must be non-chaotic in strict technical sense, 
it does exhibit "ray-chaotic footprints" in the high-frequency 
regime. These, with a few notable exceptions (see, e.g., |35|, 
|[36l ). generally result into an irregular random-like behavior. 
In such scenario, ray-based techniques are deterministic ally 
inaccurate and full-wave predictions may not be necessary 
for certain applications, while statistical models based on 
random-plane-wave (RPW) superpositions fyT\, ["381 may be 
more insightful and useful to capture and parameterize the 
essential physics. Interestingly, similar RPW models were 
independently developed within the framework of narrowband 
EM reverberation chambers |[39l . 

Against the above background results, mostly focused on 
time-harmonic excitation (for which the quantum/EM anal- 
ogy may be exploited to its fullest), this paper deals with 
short-pulsed (SP) wavepacket propagation in ray-chaotic en- 
closures. For this regime, comparatively fewer results are 
available in the quantum-physics literature (mostly dealing 
with "Loschmidt echoes" and "fidelity decay" ll40l ). and their 
translation to EM scenarios is less straightforward in view of 
the limited (to the paraxial regime) equivalence between the 
wave equation and the time-dependent Schrodinger equation. 
Among the EM results available, it is worth mentioning those 
pertaining to "one-channel" time-reversal focusing |25 | and 
related fidelity -based sensing schemes ll26l . where the ergod- 
icity properties associated with ray chaos are exploited in order 
to trade off the conventional spatial sampling with temporal 
sampling. Also worth of note is the study in iHTl . which 
illustrates the transitioning from power-law to exponential 
time-decay of the power reflected by a SP-excited one-port 
cavity with ray-chaotic geometry. 

Our interest in this subject was initially motivated by 
the speculation that pulse-reverberating (multi-echoing) ray- 
chaotic enclosures might provide spatio-temporal randomized 
field distributions of potential interest for wideband EM inter- 
ference and/or EM compatibility testbeds. Some preliminary 
studies, first via an oversimplified ray analysis PP2l and 
subsequently via rigorous numerical simulations lfT2l based 
on the finite-difference-time-domain (FDTD) method P3l . 
confirmed the viability of this intuition. In this paper, we 
develop and validate a rigorous analytical framework based 
on the extension to the SP time-domain (TD) of the RPW 
model. 
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Fig. 2. (a) Typical space-filling ray trajectory (after 250 reflections) in a 
quarter-Sinai-stadium enclosure, (b) Pseudo phase-space map displaying, for 
5000 reflections, the curvilinear abscissa [measured counterclockwise from the 
point of coordinates (x = 0,y = 0.25a), normalized to the total boundary 
length] vs. the direction cosine of the reflected ray. 



Fig. 1. Problem geometry. A 2-D (quarter) Sinai-stadium-shaped enclosure 
of characteristic size a is excited by a SP wavepacket generated by a spatio- 
temporal aperture field distribution [cf. (T)]. 



Accordingly, the rest of the paper is organized as follows. 
In Sec. ini with specific reference to a Sinai-stadium-shaped 
enclosure, we outline the physical model and related ray-based 
and wavefield observables, highlighting the peculiarities of the 
late-time response. In Sec. |III1 we develop a SP-RPW model 
which, in the parametric regime of interest, yields a spatio- 
temporal Gaussian random field whose statistics are amenable 
to analytical treatment. In Sec. |IV] we validate the SP- 
RPW model against the spatio-temporal statistical observables 
estimated from a full-wave numerical solution based on the 
FDTD method. Finally, in Sec. |V] we provide some brief 
concluding remarks and perspectives. 

II. Physical Model 

A. Geometry, Observables, and Notation 

The problem geometry is illustrated in Fig. [T] We consider 
a two-dimensional (2-D) scenario featuring a closed, hollow 
enclosure with perfectly-electric-conducting (PEC) walls and 
z— invariant quarter-Sinai-stadium geometry |44| with char- 
acteristic size a. We assume the enclosure to be excited 
by a localized SP wavepacket generated by an (electric, 
z— directed) aperture field distribution along the axis, cen- 
tered at {xA,yA), 



eAz{x,t) = A(a:: - xa)s 



t 



sin 9 a 
Co 



{x - xa) 



(1) 



with A(x) denoting a spatial window localized within the 
interval \x\ < La/2, s{t) a SP waveform localized within 
the interval \t\ < Ts/2, cq — l/^£o/^o the speed of light in 
vacuum, and 9a a possible slant angle (see Fig.[T]). This yields 
a transverse-magnetic (TM) polarized EM field distribution 

e(r,t) = e^(r, t)z, h{r,t) = hr,{r,t)x + hy{r,t)y, (2) 

with r = XX + yy, and a denoting an a— directed unit vector 
Throughout the paper, boldface symbols identify vector quanti- 
ties, lower- and upper-case symbols identify time-domain (TD) 
and frequency-domain (FD) observables, respectively, and a 



tilde (~) identifies wavenumber-domain quantities, according 
to 



fit) - 

F{u) = 

Fir) = 

F(k) = 



1 

2^ 



F{lu) exp {—iLot) duj, 
f{t) exp (iut) dt, 



^ II F(k)exp(zk.r)dk, 



F(r) exp (— ik • r) dr. 



(3a) 



(3b) 



(4a) 



(4b) 



In view of the assumed lossless character, straightforward 
application of Poynting theorem P31 yields the conservation 
of the total EM energy (per unit length in z) injected in the 
enclosure through the aperture, given by 



EM 



JJ ■WEM{r,Te)dr 
= J dt J eAix,t) X hAix,t) ■ ydx, (5) 

-oc' Va 

where Ve and Va indicate the enclosure and aperture do- 
mains, respectively, e^, the EM field at the aperture, 

WEAiir, ^) = ^ [eoe(r, t) ■ e(r, t) + fioh{r, t) ■ h(r, t)] (6) 

the EM energy density, and 

„ La sin 9a 

2 2co 

guarantees the extinction of the aperture-field transient. 



(7) 



B. Ray-Chaotic Behavior 

The quarter-Sinai-stadium geometry considered in our 
study, first introduced in B4ll . consists of a quarter stadium 
with an off-centered quarter-circle indentation (Fig. [T]), and 
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Fig. 3. SP ( normalized) waveform considered in the FDTD simulations [cf. 
127 1 and )29t . with q = 2 and ? = 0.25]. The corresponding (normalized) 
energy spectral density is shown in the inset. 



is basically a combination of the well-known Sinai (rectangle 
with central circle) [4| and Bunimovich (stadium) |5| billiard 
geometries. Similar to these canonical geometries, it exhibits 
ray-chaotic behavior but, unlike them, it is devoid of non- 
isolated periodic trajectories ("bouncing-ball" modes). 

While a comprehensive study of the ray-chaotic dynamics 
is not in order here (see, e.g., (12] for more details). Fig. 
|2] provides a basic illustration of the phenomenon. More 
specifically. Fig. |2ja) displays the space-filling properties of 
a typical ray trajectory (after 250 reflections), while Fig. 
|2b) shows a pseudo phase-space diagram, namely, for each 
reflection, the curvilinear abscissa of the impact point vs. the 
direction cosine of the reflected ray (details in the caption). 
The practically uniform covering of the phase space confirms 
the visual impression from Fig. |2|a) that (apart from zero- 
measure sets of launching conditions corresponding to iso- 
lated periodic trajectories) rays eventually approach almost 
every point arbitrarily closely and arbitrarily many times, 
with uniformly distributed arrival directions, when the number 
of reflections tends to infinity. This is markedly different 
from what observed in coordinate-separable geometries (e.g., 
rectangle, circle, ellipse), where ray trajectories tend to fill the 
enclosure in a regular fashion, with correspondingly sparsely 
populated phase spaces (see IIT2I for more details). 

Clearly, given their purely kinematical nature, the above 
considerations do not depend on the type (time-harmonic or 
SP) of excitation. 

C. Full-Wave Modeling 

For a preliminary qualitative understanding of the SP 
wavepacket propagation in the Sinai-stadium-shaped enclo- 
sure, as well as for numerical validation purposes (see Sec. 
IIVI) . we use the FDTD method ll43l . In our simulations, the 
aperture field distribution in ([T]) is implemented as a "hard 
source" li43J . with a Gaussian spatial-window function 

A(x) = exp (-^) (8) 



of effective width La ~ cl/&, and the SP waveform shown in 
Fig. |3] (see Sec. IIII-EI below for details) of effective duration 
Ts = a/(20co), chosen so as to guarantee a significant spatio- 
temporal localization on the enclosure scale (cqTs <C a). 

The computational domain is discretized with spatial steps 
Aj; = Ay = a/1200 [corresponding to ~ 38 points per center 
wavelength, and ^ 15 points per minimum wavelength, cf. the 
inset in Fig. [3, and time step At = rs/424 sa A2,/(5\/2co) 
(i.e., five-time below the Courant limit for numerical stability 
||43| ). The above spatio-temporal discretization parameters 
were chosen so as to guarantee a suitable tradeoff between 
computational affordability and accuracy (assessed via self- 
consistency tests based on energy conservation), as well as 
the need to provide a suitable number of spatio-temporal field 
samples for the statistical analysis in Sec. |IV] 

Figure |4] shows a typical temporal evolution of the wave- 
field, via four instantaneous snapshots at representative times. 
While propagating inside the cavity and undergoing mul- 
tiple reflections at the boundaries, the initial "bullet-type" 
wavepacket [Fig.|4ja)] progressively loses its shape and spatio- 
temporal localization, and its energy is eventually spread 
across the enclosure in a uniform random-like fashion [Fig. 
Ud)]. In this regime, as already noted, statistical modeling 
may provide a better problem-matched description and pa- 
rameterization, as it may effectively reproduce some coarse 
features (e.g., mean value and energy, correlation). It should 
be noted, however, that such a description may fail satisfying 
the boundary conditions, and its use should be accordingly 
limited to regions sufficiently far from the boundaries (e.g., in 
terms of the largest wavelength in the field spectrum). 

It is insightful to qualitatively compare the above be- 
havior with what occurs in "ray-regular" (i.e., non-chaotic) 
coordinate-separable geometries. To this aim. Fig. |5] shows 
the late-time snapshots pertaining to a rectangular and a semi- 
circular enclosure of same are as the quarter-Sinai-stadium 
geometry. For the rectangular case [Fig.|5ja)], a rather uniform 
and yet regular covering characterized by a checkerboard-type 
pattern is observed, which reflects the typical fixed-angle char- 
acter of the associated ray trajectories [12|. On the other hand, 
the field snapshot pertaining to the semi-circular geometry 
[Fig- 13b)] appears quite disuniform, reflecting the structure 
of typical annular-sector-filling ray trajectories bounded by 
circular caustics lfT2l . 

III. Statistical Model 

A. SP-RPW Superposition 

In ll37l . Berry conjectured that the higher-order eigenstates 
in ray-chaotic quantum billiards should be statistically in- 
distinguishable from a superposition of time-harmonic plane 
waves with random and uncorrelated directions of propagation, 
amplitudes, and phases. The RPW model, has been demon- 
strated (see, e.g., |31 1-|33 |) to effectively capture the statistical 
properties of high-order wavefunctions of strongly chaotic 
billiards in the high-frequency regime (see, e.g., 1461 . pTl 
for possible deviations due to bouncing-ball modes). In what 
follows, we extend the RPW model, perhaps for the first time, 
to the SP-TD scenario of interest here. 
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Fig. 4. FDTD-computed electric-field snapshots for the quarter-Sinai- stadium enclosure in Fig. [Hat (a) to = 3Ts, (b) to = 50Ts, (c) to = lOOTs, (d) 
to = lOOOTs. Hard-source spatio-temporal aperture field distribution as in {T), with xa = 0.75a, yA = 0.25, 6a = 30°; spatial window in ^8) with 
La = a/6; SP waveform in Fig.|3][cf 0?) and (29), with g = 2 and <; = 0.25] with duration Ts = a/(20co). 
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Fig. 5. As in Fig. |4td) (i.e., to = lOOOTs) but for geometry-separable ray- 
regular enclosures with (a) square and (b) half-circular shape, and same total 
area (not in scale). 



Let us consider a spatio-temporal domain of observation 
r=\r\<R, to<t<to + T, (9) 
with to sufficiently large so that the wavepacket energy has 



spread uniformly across the enclosure [as, e.g., in Fig. SJd)], 
and with R ^ cqTs and T ^ Tg. Moreover, we make the 
assumption that the observation points r in the domain are suf- 
ficiently far (i.e., several effective wavelengths or, equivalently, 
pulse widths cqTs away) from the enclosure boundary (see 
the discussion in Sec. Ill-CI ). We make the ansatz that the EM 
field distribution can be modeled in terms of a superposition 
of independent SP-RPWs of the form 



N 



e.ir 



t) ~ ^ a„s(i - /t„ • r - <„) , (10a) 

n=l 
1 ^ 

h(r,i) — a„s(t - K„ • r - i„) /«„ X z, (10b) 

where N, an, and t„ are independent random variables. 
More specifically, the SP-RPWs are assumed to follow a 
Poisson arrival process described by the probability function 
11481 



^T,i?.)^exp(-iV) 
iV! ' 



Pt.t?. {N) 



(11a) 



with Nt,r representing the average number of arrivals in 
the spatio-temporal domain of observation (|9]l. The amplitude 
coefficients a„ are uncorrected, zero-mean, unit- variance ran- 
dom variables. 



(lib) 



with (•) denoting the statistical (ensemble) expectation value, 
5nm the Kronecker delta, and no specific assumption made on 
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higher-order correlations. The "normahzed" wavevectors are 
expressed by 

«ri = Cp"^ (cos0„x + sin(/)„y), (11c) 

where {</)„ } is a set of uniformly, independently distributed 
random variables within the interval [—tt, tt]. The time-delays 
tn are uniformly, independently distributed within the interval 
[to — R/co, to + T+R/co], which guarantees causal connection 
with the spatio-temporal observation domain in (|9]l. 

B. Background Results on Poisson Impulse Noise Processes 

In order to understand the statistical properties of the SP- 
RPW superpositions in ( fTOl i. we capitalize on a body of well- 
established technical results on Poisson impulse noise pro- 
cesses. In particular, following I.48J, a critical nondimensional 
parameter for establishing the statistical structure is 

Nt,rTs 



with 



7 ^ 



2R 

Co 



(12) 



which basically represents the average number of SP-RPWs 
per unit time at a fixed observation point. Intuitively, small 
values of 7 (up to ~ 10 I.48J) imply weakly overlapping 
SP-RPWs and, consequently, a strong dependence of the 
statistics on the SP waveform s{t) as well as on the specific 
distributions of the random parameters a„, k„ and t„, with a 
high probability associated with zero or small amplitude values 
||48|. This regime may somehow resemble the intermediate- 
time response of the enclosure, when the wavepacket energy 
is only partially spread across the enclosure [see, e.g.. Fig. 
|4tc)], and is not of direct interest for this investigation. In fact, 
the late-time response [cf. Fig. SJd)] of interest here is more 
intuitively associated with strongly overlapping SP-RPWs, and 
hence large values of 7. In this regime, as a consequence of 
the central limit theorem |48J, the statistics approach Gaussian 
behavior More specifically ll48l . nearly-Gaussian behaviors 
(with correction terms of order 7"^/^ or 7"^) should be 
expected for 7 values ranging from ^ 10 to ^ 10'', whereas 
7 > 10^ should ensure Gaussian behavior irrespective of the 
SP waveform shape and specific parameter distributions, with 
the parameter 7 appearing as a scale factor in the probability 
density function (PDF). In this regime, the (insofar unknown) 
parameter Nt.r may be related to the actual physical scenario 
via simple energy considerations (see Sec. IIII-DI below). 

C. Spatio-Temporal Gaussian Random Field 

It follows from the considerations above that the late-time 
SP response of our ray-chaotic enclosure may be effectively 
modeled by a spatio-temporal Gaussian random field |49l. 
Accordingly, given a JC-dimensional array of spatio-temporal 
field samples 

e^=[ez{ri,Ti) e^{T2,T2) ■■■ e^{rK,TK)f , (13) 
with denoting the transpose, the associated joint PDF is 



given by [|49i 



{2tt)K det(E) 



■ exp 



(14a) 



a ^ d^) ' 1 ^ ( (££ " ■ (££ " /f) 



(14b) 



representing the expectation-value vector and covariance ma- 
trix, respectively. It can readily be shown (see Appendix A for 
details) that 

(e.(r,0)=0, (15) 

and therefore the Gaussian joint PDF in ( I14al ) is completely 
determined by the knowledge of the spatio-temporal correla- 
tion 

c(p, r) = {e, (r, t) e, {r + p,t + t)) , (16) 

where stationarity in space and time (i.e., independence on 
r and t) is anticipated. Alternatively, the correlation in ( fT6b 
may also be calculated via a spectral (Wiener-Khinchin-type) 
representation 1491 , 



c(p, r) = ^-g J duj JJ dk CeiK w) exp [i (k • p - iujr)] , 
—00 —00 

(17) 

where Ce(k, w) is the power spectral density (PSD). Simple 
heuristic arguments would suggest a generalization of the 
known FD results ifJTl in terms of a separable PSD, 



Co 



(18) 



where S is a multiplicative constant, the Dirac-delta spectral 
kernel enforces the plane-wave dispersion relation in vacuum, 
and [^(a;)^ accounts for the energy spectral density of the 
SP waveform s{t). In view of the assumed structure, the 
wavevector integral in ([TtI i can be performed analytically 
via a simple change of variables to cylindrical coordinates, 
k = wcq "'"(cos^sx + sin^sy), p = /^(cos (^qx + sin (^oy), 
yielding 



c{p,t) - 



dLL>\S{LL>)\ exp(— iwr) 



exp 



I cos [ip — (po) 

Co 



dip 



47r2 



I^HlVo 



Co 



exp (— iwr) dw, (19) 



where Jq denotes a zeroth-order Bessel function lISOl . We note 
that the result in ( fT9l ) does not depend on the direction of 
the spatial displacement vector p, thereby implying isotropy 
of the spatio-temporal Gaussian random field. Moreover, Eq. 
(T% reproduces the well-known Jo-type FD form |37| in the 
time-harmonic limit. 

As shown in Appendix B, the above heuristic-based result 
is consistent with first-principle calculations [i.e., with (fTSI l 
and (fTotl. In particular, as also shown in Appendix B, we 
can relate the multiplicative constant S in iT% to the critical 
parameters in the SP-RWP model, viz.. 



(20) 
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whereby the anticipated appearance of the critical parameter 
7 in (fTZt as a scaling factor in the Gaussian PDF. 

D. Link with Underlying Physics 

A link between the SP-RPW model and the underlying 
physics can be established via simple energy considerations. 
First, via straightforward vector algebra, we note that, in spite 
of the random character of the normalized wavevector /«„ in 
([Tol l, the quantity 

(k„ X z) • {Kn X z) = Cq ^ (21) 

is always deterministic. This implies that the spatio-temporal 
correlation for the magnetic field is simply related to that of 
the electric field in (fTST i via 



(h (r,i) •h(r + p,i + T)) 



(22) 



where 770 = yyUo/so is the vacuum characteristic impedance. 
Moreover, Eq. (|22] | also establishes a simple connection with 
the expectation value of the EM energy density in (|6]l. 



{wEM{r,t)) = eoc(0,0), 



(23) 



which, in view of spatio-temporal stationarity, is in turn simply 
related to the total EM energy (per unit length) in (|5]l via 



Ae 



(24) 



where Ae is the area of the enclosure. Accordingly, we can 
rewrite the spatio-temporal correlation in ( fT9l ) as 



c(p,t) 



where 



/ £em 
XSoAe 



15(^)1^0 ( ^ ) exp(-iwT)dw, (25) 



\SiL0)\' 



(26) 



represents the normalized energy spectral density of the SP 
waveform. The expression in ( |25] l relates the statistical prop- 
erties of the SP-RPW model to simple geometrical and energy 
invariants of the physical model. 

E. Analytically Workable Examples 

It is instructive to work out a class of examples for which 
the integral in (l25T l can be calculated analytically in closed 
form. To this aim, we consider a class of SP waveforms 

* 2 

where is a scaling parameter (chosen so as to ensure that 
the effective waveform duration is Tg within a suitably small 
error), the superfix ^''^ indicates the th derivative, and 



Sit) 



s(t) = Re 



(27) 




lm{t) = 



(28) 



is the analytic delta function pTl (with V denoting the 
principal value). It can readily be shown that 



5H-(-zc.)«cxp(-|u;|^ 



(29) 



which, substituted in yields (cf. Eq. (2.12.8.4) in [52)) 

'(2<Z-I-1) 



^^^'"^ = soAe 



X P, 



(30) 



where denotes the mth-degree Legendre polynomial (cf. 
Sec. 8 in ISl), and 

e= f rS^- . (31) 

Co (<jTs + it) 

is a complex, nondimensional spatio-temporal variable. Partic- 
ularly simple and insightful is the special case q = (analytic 
delta), for which 

1-1/2^ 



_ I ^Ts£em\ 
sqAe J 



Re 



it) 



whose 1-D cuts at r = and p = Q can be written as 

Pc^EM 



c(p,0) 



eoAE^/il? 



Pi 



c(0,t) = 



-EM 



(32a) 



(32b) 



(32c) 



EoAe (r2 + T2) ' 

with pc = \/icQ<^Ts and Tc = <,Ts representing the correlation 
length and time (3dB below), respectively. This simple exam- 
ple completes the illustration of how the statistical properties 
of the late-time field distribution, modeled as a spatio-temporal 
Gaussian random field, are related to (and may be potentially 
controlled via) a minimal number of generic geometrical and 
physical parameters. 

IV. Numerical Validation 

A. Generalities and Observables 

We now move on to validating the SP-RPW model devel- 
oped in Sec. against the statistical observables that can be 
directly estimated from the FDTD-computed late-time wave- 
field distributions [cf. Fig. IHd)]. In order to facilitate direct 
comparison with the theoretical predictions, we introduce a 
suitably normalized field, 



eoAi 



8 



EM 



(33) 



for which the predicted statistics are standard (i.e., zero-mean, 
unit-variance) Gaussians. In ( |33] ), Ae = (3 -f 87r)a2/16 (cf. 
Fig. [T]), the wavelfield ez is computed as detailed in Sec. 
III-CI and the EM energy (per unit length) £em is computed 
numerically via spatial integration of the EM energy density 
[cf. (|5])-(|7|]. Note that, in order to avoid possible direct-current 
field offset issues in the FDTD simulations |43|, we prefer to 
work with the zero-mean SP waveform in Fig. [51 rather than 
the simplest analytic-delta (q — 0) case in ( |32] ). 
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Fig. 6. Gaussian probability plots (blue-solid curves) for local (a) temporal 
and (b) spatial statistics of the normalized field in B3t estimated from a 
single FDTD-computed realization. Local temporal statistics are estimated 
observing the wavefield at a fixed point (Xc = 0.602a, j/c = 0.121a) over 
an interval of duration 20Ts (i.e., M = 8478 samples) starting at to = 
lOOOTs (estimated mean and variance are —0.0035 and 1.07, respectively). 
Local spatial statistics are estimated over a AcqTs x AcqTs square domain 
(i.e., M = 57121 samples) centered at point (Xc,yc) at time to = lOOOTs 
(estimated mean and variance are —0.0236 and 1.057, respectively). As a 
reference, also shown (red-dashed straight line) is the behavior of a standard 
(zero-mean, unit-variance) Gaussian distribution. 

Assuming ergodicity, the estimation of the statistical (en- 
semble) observables of interest is based on single spatial or 
temporal wavefield realizations. Moreover, in order to reveal 
possible variations (in time, and across the enclosure), we 
are interested in local statistics, and accordingly consider 
observation domains with sufficiently small size (compared 
with the overall spatio-temporal observation domain), and yet 
sufficiently large (on the SP wavepacket scale) so as to guar- 
antee meaningful statistics. More specifically, unless otherwise 
stated, local temporal statistics are computed over time-series 
of size 20Ts (i.e., 8478 samples) at a fixed observation point, 
while local spatial statistics are computed over square domains 
of size 4:CoTs x AcqTs (i.e., 57121 samples) around a given 
point. 

B. Representative Results 

For a first qualitative analysis, we show in Fig. |6] the Gaus- 
sian probability plots pertaining to point (xc = 0.602a, j/c = 
0.121a) (magenta down-triangle marker in the inset of Fig. 
|8]l at a late time to = lOOOTs, when the SP wavepacket 
energy has spread uniformly across the enclosure [cf. Fig. 
|4jd)]. More specifically, with reference to the local temporal 
[Fig-ISa)] and spatial [Fig.|6jb)] distributions, the plots display 
the cumulative counts pertaining to the normalized field in (l33T l 
using a Gaussian probability scale, so that standard-Gaussian- 
distributed data sets would yield a reference straight line 
(shown as red-dashed). This facilitates direct visual assessment 
of the Gaussianity of the data, or otherwise the presence of 
short/long tails and/or left/right skewness. As evident from 
Fig. |6] at the chosen late time, both the temporal and spa- 
tial wavefield distributions follow quite closely the reference 
standard-Gaussian behavior predicted by the SP-RPW model. 



Fig. 7. As in Fig. |6] but for earlier-time (to = lOOTs) observations, (a) 
Local temporal statistics (estimated mean and variance are 0.0018 and 3.389, 
respectively), (b) Local spatial statistics (estimated mean and variance are 
-0.00649 and 0.448, respectively). 
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Fig. 8. Normalized KS distance in (34} pertaining to temporal field 
observations at five randomly-chosen points (shown as different markers in 
the insets) over an interval of duration 20Ts (i.e., M = 8478 samples) for 
different starting times to- 

Figure |7] shows the Gaussian probability plots pertaining 
to an earlier-time [t — lOOTs, cf. Fig. SJc)] response at the 
same observation point, which turn out to be markedly non- 
Gaussian (see also the estimated mean and variance values 
given in the caption). Diverse types of distributions, not shown 
here for brevity, were observed for different observation points, 
confirming the general trend of markedly non-Gaussian behav- 
ior with significant parameter dependence on the observation 
point. 

For a more quantitative assessment and compact visual- 
ization of the temporal and spatial fluctuations of the local 
statistics, it is expedient to introduce a measure of distribution 
distance, namely, the Kolmogorov-Smirnov (KS) |,54J , 

Dm^&M^m{X)-HX)\, (34) 

X 

where <I>m is the cumulative density function (CDF) numeri- 
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Fig. 9. As in Fig. [8] but for spatial observations over 4coTs X 4coTs square 
domains (i.e., M = 57121 samples) centered at the same randomly-chosen 
points. 



cally estimated over M samples, and 



(35) 



is the standard Gaussian CDF. Figure |8] shows the evolution 
of the (normalized) KS distance in ( |34] | pertaining to the 
local temporal statistics at five randomly-chosen observation 
points (shown as different markers in the inset) with distances 
> 2cqTs from the enclosure boundary, as a function of the 
observation time. It is evident that, starting from rather large 
and strongly position-dependent values at early times, the 
KS distance progressively decreases with increasing obser- 
vation times, converging to much smaller (about an order 
of magnitude) and uniform values for sufficiently late times. 
For quantitative reference, we recall that, in typical KS- 
based goodness-of-fit tests, the null hypothesis (i.e., Gaussian 
distribution) is rejected with a 0.01 significance level (i.e., 99% 
confidence level) if ll54l 

^/MDM > 1.63. (36) 

Therefore, although a strict goodness-of-fit assessment is not 
in order here, we can safely conclude that the SP-RPW- 
based Gaussian prediction seems to model reasonably well, 
and fairly uniformly across the enclosure, the late-time local 
temporal statistics. Similar observations and conclusions hold 
for the local spatial statistics over neighborhoods of the same 
observation points, as shown in Fig. |9] 

Figures [T0la)- [T0f c) show three representative snapshots of 
the late-time spatio-temporal correlation [estimated within a 
neighborhood of point (xc — 0.602a, t/c — 0.121a)] for 
different values of the time lag r, with the corresponding 
1-D cuts compared in Figs. [TOld)- [TO[ f) with the theoretical 
predictions from the SP-RPW model [cf. with q = 2 
and ? = 0.25]. Apart from weak fluctuations attributable to 
edge effects induced by the limited extent (Ac^Ts x AcqTs) 
of the observation domain, the dominant features appear 
rather isotropic, and well captured by the SP-RPW prediction. 



Finally, Fig. [TTI shows the late-time temporal (i.e., zero spatial 
displacement) correlation estimated at a fixed observation 
point over an interval of duration 200Ts, compared with the 
SP-RPW prediction. Once again, a quite good agreement 
is observed. Qualitatively similar results were consistently 
observed for different observation points across the enclosure. 

The above results clearly demonstrate that the SP-RPW 
model developed in Sec.|lII]effectively captures and parameter- 
izes the essential statistical features of the irregular random- 
like regime that characterizes the late-time response of our 
ray-chaotic enclosure in the SP regime. 

It is interesting to compare our SP-TD results with the 
well-known time-harmonic case. In that case, the local spatial 
statistics (estimated over regions of size much larger than 
the wavelength, but much smaller than the enclosure size) 
of high-order eigenfields approach Gaussian behavior, with 
spatial correlation approaching a Bessel function (Jq) form, 
which is well captured by time-harmonic RPW superpositions 
ll37l . In our SP-TD case, the late-time local spatio-temporal 
statistics (estimated over regions of size much larger than 
the SP scale, but much smaller than the observation domain) 
approach Gaussian behavior, with spatio-temporal correlation 
depending on the SP waveform and coarse (geometrical and 
energy) parameters via 



V. Conclusions and Perspectives 

In this paper, we have developed and validated a SP- 
RPW theoretical framework for the statistical modeling of the 
SP wavepacket propagation in a ray-chaotic enclosure. With 
specific reference to the quarter-Sinai-stadium geometry, our 
results confirm that the late-time response of the enclosure can 
be effectively modeled as a spatio-temporal random Gaussian 
field, whose parameters are simply related to coarse geomet- 
rical quantities, energy invariants, and the SP waveform. 

The peculiar (and, to a certain extent, controllable) features 
of such wavefield distribution may find important applications 
to wideband EM interference and/or EM compatibility testbeds 
(whereby an equipment under test may be subjected to a "pulse 
shower" illumination with characteristics largely independent 
of its location, orientation, shape and constitutive properties), 
as well as the synthetic emulation of multi-path wireless 
channels. Within this framework, our results, based on the 
idealized assumption of the absence of losses, may be extended 
to a more realistic low-loss scenario by assuming adiabatic 
variations of the statistics as a consequence of the slow (on 
the SP time-scale) energy decay. Also of great interest is 
the exploration of electronic-stirring mechanisms based on 
narrowband-excited reverberating chambers with ray-chaotic 
geometries. 
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Fig. 10. Top panels: (a), (b), (c) Spatio-temporal correlation snapshots (as a function of p = p^x + pyy) at time to = lOOOTs, estimated over a 
icoTs X AcoTs square domain (i.e., M = 57121 samples) centered at point (Xc = 0.602a, yc = 0.121a), for r = 0, O.25T3 , 0.5Ts, respectively. Bottom 
panels: (d), (e), (f) Corresponding 1-D cuts (blue-solid) aX py = compared with SP-RPW prediction (red-dashed) from 130) (with q = 2 and <; = 0.25). 
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Fig. 11. Temporal (p = 0) coirelation (blue-solid) estimated from a 
wavefield observation at a fixed point (xc = 0.602a, j/c = 0.121a) over 
an interval of duration 200Ts (i.e., AI = 84780 samples) starting at 
to = lOOOTs, compared with SP-RPW prediction (red-dashed) from (30) 
(with g = 2 and f = 0.25). 
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Appendix A 
Pertaining to Eq. (fTsI) 

In order to calculate the electric-field expectation value 



N 



{e^ (r, t)) ^ ans{t - /«„ • r - t„) 



(37) 



\n=l 



we recall Wald's identity |55| on the expectation value of 
the sum of a random number N of finite-mean, identically 
distributed random variables X„ (independent of N), viz.. 



N 



Y,Xn) = {N) (X„) 



(38) 



Straightforward application of ( |38] | to ( l37b . and simple fac- 
torization (in view of the assumed statistical independence), 
yields 

(e, (r, t)) - Nt,r (a„) {s{t - k„ • r - t„)) , (39) 
from which Eq. dTsl l follows immediately recalling jllbl) . 

Appendix B 
Pertaining to Eq. (|201) 



The first-principle calculation of the spatio-temporal cor- 
relation via ( fT6b with ( fTOl i is most easily carried out in the 
FD, using the spectral representation [via (O] of the SP-RPW 
model 

N °° 
1 , ^ f 

62 (r, i) — ^ a„ / S{uj)exp[iuj {Kn ■ r + tn - t)]duj. 
^ "=i -00 

(40) 

It is also expedient to rewrite ( fTSI ) as 

c(p,t) = {e,{v,t)el{v + p,t + T)), (41) 
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where the "dummy" (in view of the real-valued wavefield) jllcl l. whence Il50ll 
complex conjugation simplifies the following technical devel- n 

opments. By substituting ( |40] i in ( 1411 1. and rearranging terms, u^„tA, ,^ „\\ _ ^ f ^^-r^ iojp cos {(j)n — ffio) 

[eXp[lUJhin ■ Pj) — r, I t-XP 

we obtain ^tt J 



Co 



N 



\rum=l _^ 



Finally, recalling ( fTSI ). we obtain 

X exp [iwi (/«„ • r + t„)] 

oo 

V c(p,r)~ — — / \S[uj)\ Jq\ — exp (-iwr) dw, 



^4:TtTsJ J V Co 

X exp [-1 (wi - W2) t + ia;2r] ^ . (42) -00 ^^^^ 

which is structurally identical, apart from the multiplicative 
Factorizing out (in view of the assumed statistical indepen- constant, to the result in (HU) obtained via a heuristic argument, 
dence) the correlation term < a„a,„ > and recalling ([TTbl), From (O and (gT]), Eq. ^ readily follows. 
Eq. dUli yields 

/ AT °° 

'^('°''^) 4^ ( X] // dujiduj2S{uji)S*{uj2) 



— 00 



X exp [i (wi - UJ2) tn] 

X exp [-i (w2 - wi) (t - K„ • r)] 

X exp [iuj2 {t- Kn- p)] ) ■ (43) 



Applying once again Wald's identity ( I38I I. and recalling that 

(exp[i (wi -W2)t„]) 
1 



exp [z (wi - a;2) in] dtn 
exp [i (wi -W2)to] 



exp [i {UJI — UJ2) tn] dtn 



7rexp[i (wi ~W2)to] x .-.s 

d[uji~'i^2), (44) 









Co J 



with the last asymptotic equality valid for R ^ cqTs and 
T ^ Ts (so as to approximately extend the integration limits 
to Too), we obtain 

, , Nt.r 

C(p, T 



4.(T+f) 

00 

X j \S{ijj)f' (exp {iuJKn ■ p)) exp {—iajr) dw. (45) 
—00 

The expectation value with respect to k„ can be easily 
computed by letting p = p{cosLpoX + s'mipoy) and recalHng 



